﻿using System;
using System.Collections.Generic;
using System.Drawing;
using System.IO;
using System.Text;

namespace WinForm
{
    public class TIN
    {
        private List<TIN_Point> points = new List<TIN_Point>();//存储离散点数据
        private List<TIN_Line> lines = new List<TIN_Line>();//存储基线数据
        private List<TIN_Tri> tris = new List<TIN_Tri>();//存储三角形数据
        private List<Contour> contours = new List<Contour>(); //存储等高线数据

        //获取离散点数据
        public bool input_points(string data)
        {
            string[] info = data.Split(',');
            try
            {
                TIN_Point item = new TIN_Point();
                item.X = Convert.ToDouble(info[0]);
                item.Y = Convert.ToDouble(info[1]);
                item.Z = Convert.ToDouble(info[2]);
                points.Add(item);
                return true;
            }
            catch
            {
                return false;
            }
        }
        
        //根据最大角原则获取第三点
        private int Max(TIN_Line line, List<int> usablepoints)
        {
            double max = TOOL.Angle(points[line.P0], points[line.P1], points[usablepoints[0]]);
            int index = usablepoints[0];
            for (int i = 1; i < usablepoints.Count; i++)
            {
                double temp = TOOL.Angle(points[line.P0], points[line.P1], points[usablepoints[i]]);
                if (max <= temp)
                {
                    max = temp;
                    index = usablepoints[i];
                }
            }
            return index;
        }

        //生成三角形邻接关系
        private void Topology()
        {
            for (int i = 0; i < tris.Count; i++)
            {
                //循环判断Tris[i]的邻接三角形
                for (int j = 0; j < tris.Count; j++)
                {
                    if (i != j)
                    {
                        if (TOOL.Equaltion2(tris[i], tris[j]))
                        {
                            int sign = TOOL.Coincide(tris[i], tris[j]);
                            tris[i].Adjoin[sign] = j;
                        }
                    }
                }
                //判断是否为边界三角形
                if (Array.IndexOf(tris[i].Adjoin, -1) != -1)
                {
                    tris[i].Edge = true;
                }
                else tris[i].Edge = false;
            }
        }

        //创建不规则三角网
        public void create_tin()
        {
            if (points.Count < 3)
            {
                return;
            }
            # region 生成第0个三角形

            /*--生成第0个三角形--*/
            TIN_Tri tri = new TIN_Tri();
            /*寻找与第0点距离最近点*/
            double mindis = TOOL.Distance(points[0], points[1]);
            int index = 1;
            for (int i = 1; i < points.Count; i++)
            {
                double temp = TOOL.Distance(points[0], points[i]);
                if (temp < mindis)
                {
                    mindis = temp;
                    index = i;
                }
            }

            /*初始基线L0*/
            TIN_Line line0 = new TIN_Line();
            line0.P0 = 0;
            line0.P1 = index;
            line0.Count = 1;
            lines.Add(line0);
            tri.L0 = lines.Count - 1;

            /*生成可用点集*/
            List<int> usablepoints = new List<int>();
            for (int i = 1; i < points.Count; i++)
            {
                if (i != index)
                {
                    usablepoints.Add(i);
                }
            }

            /*根据张角最大原则获取第三点*/
            index = Max(line0, usablepoints);

            /*生成基线L1*/
            TIN_Line line1 = new TIN_Line();
            line1.P0 = line0.P0;
            line1.P1 = index;
            line1.Count = 1;
            lines.Add(line1);
            tri.L1 = lines.Count - 1;

            /*生成基线L2*/
            TIN_Line line2 = new TIN_Line();
            line2.P0 = line0.P1;
            line2.P1 = index;
            line2.Count = 1;
            lines.Add(line2);
            tri.L2 = lines.Count - 1;
            //添加第0个三角形至链表中
            tris.Add(tri);
            #endregion 
            
            /*对于第K个三角形为拓展三角形时*/
            int k = 0;
            while (true)
            {
                #region 处理第K个三角形的L0边
                /*处理第K个三角形的L0边*/
                if (lines[tris[k].L0].Count < 2)
                {
                    /*获取第K个三角形，不是L0边上的端点*/
                    int test_point = -1;
                    if (lines[tris[k].L0].P0 == lines[tris[k].L1].P0 ||
                        lines[tris[k].L0].P1 == lines[tris[k].L1].P0)
                    {
                        test_point = lines[tris[k].L1].P1;
                    }
                    if (lines[tris[k].L0].P0 == lines[tris[k].L1].P1 ||
                        lines[tris[k].L0].P1 == lines[tris[k].L1].P1)
                    {
                        test_point = lines[tris[k].L1].P0;
                    }

                    /*获取该端点在竖直方向到L0线段的距离*/
                    double sign1 = TOOL.Sign(points[lines[tris[k].L0].P0], points[lines[tris[k].L0].P1], points[test_point]);

                    /*清空可用点集*/
                    usablepoints.Clear();

                    /*获取可用点集*/
                    for (int i = 0; i < points.Count; i++)
                    {
                        if (lines[tris[k].L0].P0 != i &&
                            lines[tris[k].L0].P1 != i &&
                            lines[tris[k].L1].P0 != i &&
                            lines[tris[k].L1].P1 != i &&
                            lines[tris[k].L2].P0 != i &&
                            lines[tris[k].L2].P1 != i)
                        {
                            double sign2 = TOOL.Sign(points[lines[tris[k].L0].P0], points[lines[tris[k].L0].P1], points[i]);
                            if (sign1 * sign2 < 0)
                            {
                                usablepoints.Add(i);
                            }
                        }
                    }

                    /*生成下一个三角形*/
                    if (usablepoints.Count > 0)
                    {
                        /*根据张角最大原则获取第三点*/
                        index = Max(lines[tris[k].L0], usablepoints);

                        /*生成临时L1、L2基线*/
                        line1 = new TIN_Line();
                        line1.P0 = lines[tris[k].L0].P0;
                        line1.P1 = index;
                        line2 = new TIN_Line();
                        line2.P0 = lines[tris[k].L0].P1;
                        line2.P1 = index;

                        /*标记临时基线是否存在的标志*/
                        int lineIdx1 = -1;
                        int lineIdx2 = -1;
                        int count1 = 0;
                        int count2 = 0;

                        /*判断基线数据集中是否存在这两条边*/
                        for (int i = 0; i < lines.Count; i++)
                        {
                            if (TOOL.Equaltion(line1, lines[i]))
                            {
                                lineIdx1 = i;
                                count1 = lines[i].Count;
                            }
                            else if (TOOL.Equaltion(line2, lines[i]))
                            {
                                lineIdx2 = i;
                                count2 = lines[i].Count;
                            }
                        }

                        /*生成最终L1、L2基线*/
                        if (count1 < 2 && count2 < 2)
                        {
                            /*新三角形的L0基线*/
                            tri = new TIN_Tri();
                            tri.L0 = tris[k].L0;

                            /*新三角形的L1基线*/
                            if (count1 == 0)
                            {
                                line1.Count = 1;
                                lines.Add(line1);
                                tri.L1 = lines.Count - 1;
                            }
                            else
                            {
                                lines[lineIdx1].Count = 2;
                                tri.L1 = lineIdx1;
                            }

                            /*新三角形的L2基线*/
                            if (count2 == 0)
                            {
                                line2.Count = 1;
                                lines.Add(line2);
                                tri.L2 = lines.Count - 1;
                            }
                            else
                            {
                                lines[lineIdx2].Count = 2;
                                tri.L2 = lineIdx2;
                            }

                            tris.Add(tri);
                            /*第K个三角形的L0基线拓展三角形生成结束*/
                        }
                    }
                    /*不论是否生成新三角形，最终都要将 第K个三角形的L0边 使用次数改为2*/
                    lines[tris[k].L0].Count = 2;
                }
                #endregion

                #region 处理第K个三角形的L1边
                /*处理第K个三角形的L1边*/
                if (lines[tris[k].L1].Count < 2)
                {
                    /*获取第K个三角形，不是L1边上的端点*/
                    int test_point = -1;
                    if (lines[tris[k].L1].P0 == lines[tris[k].L0].P0 ||
                        lines[tris[k].L1].P1 == lines[tris[k].L0].P0)
                    {
                        test_point = lines[tris[k].L0].P1;
                    }
                    if (lines[tris[k].L1].P0 == lines[tris[k].L0].P1 ||
                        lines[tris[k].L1].P1 == lines[tris[k].L0].P1)
                    {
                        test_point = lines[tris[k].L0].P0;
                    }

                    /*获取该端点在竖直方向到L1线段的距离*/
                    double sign1 = TOOL.Sign(points[lines[tris[k].L1].P0], points[lines[tris[k].L1].P1], points[test_point]);

                    /*清空可用点集*/
                    usablepoints.Clear();

                    /*获取可用点集*/
                    for (int i = 0; i < points.Count; i++)
                    {
                        if (lines[tris[k].L0].P0 != i &&
                            lines[tris[k].L0].P1 != i &&
                            lines[tris[k].L1].P0 != i &&
                            lines[tris[k].L1].P1 != i &&
                            lines[tris[k].L2].P0 != i &&
                            lines[tris[k].L2].P1 != i)
                        {
                            double sign2 = TOOL.Sign(points[lines[tris[k].L1].P0], points[lines[tris[k].L1].P1], points[i]);
                            if (sign1 * sign2 < 0)
                            {
                                usablepoints.Add(i);
                            }
                        }
                    }

                    /*生成下一个三角形*/
                    if (usablepoints.Count > 0)
                    {
                        /*根据张角最大原则获取第三点*/
                        index = Max(lines[tris[k].L1], usablepoints);

                        /*生成临时L1、L2基线*/
                        line1 = new TIN_Line();
                        line1.P0 = lines[tris[k].L1].P0;
                        line1.P1 = index;
                        line2 = new TIN_Line();
                        line2.P0 = lines[tris[k].L1].P1;
                        line2.P1 = index;

                        /*标记临时基线是否存在的标志*/
                        int lineIdx1 = -1;
                        int lineIdx2 = -1;
                        int count1 = 0;
                        int count2 = 0;

                        /*判断基线数据集中是否存在这两条边*/
                        for (int i = 0; i < lines.Count; i++)
                        {
                            if (TOOL.Equaltion(line1, lines[i]))
                            {
                                lineIdx1 = i;
                                count1 = lines[i].Count;
                            }
                            else if (TOOL.Equaltion(line2, lines[i]))
                            {
                                lineIdx2 = i;
                                count2 = lines[i].Count;
                            }
                        }

                        /*生成最终L1、L2基线*/
                        if (count1 < 2 && count2 < 2)
                        {
                            /*新三角形的L0基线*/
                            tri = new TIN_Tri();
                            tri.L0 = tris[k].L1;

                            /*新三角形的L1基线*/
                            if (count1 == 0)
                            {
                                line1.Count = 1;
                                lines.Add(line1);
                                tri.L1 = lines.Count - 1;
                            }
                            else
                            {
                                lines[lineIdx1].Count = 2;
                                tri.L1 = lineIdx1;
                            }

                            /*新三角形的L2基线*/
                            if (count2 == 0)
                            {
                                line2.Count = 1;
                                lines.Add(line2);
                                tri.L2 = lines.Count - 1;
                            }
                            else
                            {
                                lines[lineIdx2].Count = 2;
                                tri.L2 = lineIdx2;
                            }

                            tris.Add(tri);
                            /*第K个三角形的L1基线拓展三角形生成结束*/
                        }
                    }
                    /*不论是否生成新三角形，最终都要将 第K个三角形的L1边 使用次数改为2*/
                    lines[tris[k].L1].Count = 2;
                }
                #endregion

                #region 处理第K个三角形的L2边
                /*处理第K个三角形的L2边*/
                if (lines[tris[k].L2].Count < 2)
                {
                    /*获取第K个三角形，不是L2边上的端点*/
                    int test_point = -1;
                    if (lines[tris[k].L2].P0 == lines[tris[k].L0].P0 ||
                        lines[tris[k].L2].P1 == lines[tris[k].L0].P0)
                    {
                        test_point = lines[tris[k].L0].P1;
                    }
                    if (lines[tris[k].L2].P0 == lines[tris[k].L0].P1 ||
                        lines[tris[k].L2].P1 == lines[tris[k].L0].P1)
                    {
                        test_point = lines[tris[k].L0].P0;
                    }

                    /*获取该端点在竖直方向到L2线段的距离*/
                    double sign1 = TOOL.Sign(points[lines[tris[k].L2].P0], points[lines[tris[k].L2].P1], points[test_point]);

                    /*清空可用点集*/
                    usablepoints.Clear();

                    /*获取可用点集*/
                    for (int i = 0; i < points.Count; i++)
                    {
                        if (lines[tris[k].L0].P0 != i &&
                            lines[tris[k].L0].P1 != i &&
                            lines[tris[k].L1].P0 != i &&
                            lines[tris[k].L1].P1 != i &&
                            lines[tris[k].L2].P0 != i &&
                            lines[tris[k].L2].P1 != i)
                        {
                            double sign2 = TOOL.Sign(points[lines[tris[k].L2].P0], points[lines[tris[k].L2].P1], points[i]);
                            if (sign1 * sign2 < 0)
                            {
                                usablepoints.Add(i);
                            }
                        }
                    }

                    /*生成下一个三角形*/
                    if (usablepoints.Count > 0)
                    {
                        /*根据张角最大原则获取第三点*/
                        index = Max(lines[tris[k].L2], usablepoints);

                        /*生成临时L1、L2基线*/
                        line1 = new TIN_Line();
                        line1.P0 = lines[tris[k].L2].P0;
                        line1.P1 = index;
                        line2 = new TIN_Line();
                        line2.P0 = lines[tris[k].L2].P1;
                        line2.P1 = index;

                        /*标记临时基线是否存在的标志*/
                        int lineIdx1 = -1;
                        int lineIdx2 = -1;
                        int count1 = 0;
                        int count2 = 0;

                        /*判断基线数据集中是否存在这两条边*/
                        for (int i = 0; i < lines.Count; i++)
                        {
                            if (TOOL.Equaltion(line1, lines[i]))
                            {
                                lineIdx1 = i;
                                count1 = lines[i].Count;
                            }
                            else if (TOOL.Equaltion(line2, lines[i]))
                            {
                                lineIdx2 = i;
                                count2 = lines[i].Count;
                            }
                        }

                        /*生成最终L1、L2基线*/
                        if (count1 < 2 && count2 < 2)
                        {
                            /*新三角形的L0基线*/
                            tri = new TIN_Tri();
                            tri.L0 = tris[k].L2;

                            /*新三角形的L1基线*/
                            if (count1 == 0)
                            {
                                line1.Count = 1;
                                lines.Add(line1);
                                tri.L1 = lines.Count - 1;
                            }
                            else
                            {
                                lines[lineIdx1].Count = 2;
                                tri.L1 = lineIdx1;
                            }

                            /*新三角形的L2基线*/
                            if (count2 == 0)
                            {
                                line2.Count = 1;
                                lines.Add(line2);
                                tri.L2 = lines.Count - 1;
                            }
                            else
                            {
                                lines[lineIdx2].Count = 2;
                                tri.L2 = lineIdx2;
                            }

                            tris.Add(tri);
                            /*第K个三角形的L2基线拓展三角形生成结束*/
                        }
                    }
                    /*不论是否生成新三角形，最终都要将 第K个三角形的L2边 使用次数改为2*/
                    lines[tris[k].L2].Count = 2;
                }
                #endregion

                k++;
                if (k == tris.Count)
                {
                    break;
                }
            }

            //生成三角形邻接表
            Topology();
        }

        //找到三角形中另一个等值点所在边的编号和位置(位置：该边在L0\L1\L2的位置，有0、1、2三种可能，方便找邻接三角形)
        private int lineindex(int edgeline,TIN_Tri A,double h,ref int m)
        {
            if (edgeline != A.L0)
            {
                TIN_Point P0 = points[lines[A.L0].P0];
                TIN_Point P1 = points[lines[A.L0].P1];
                if (TOOL.Existence(P0.Z, P1.Z, h))
                {
                    m = 0;
                    return A.L0;
                }
            }
            if (edgeline != A.L1)
            {
                TIN_Point P0 = points[lines[A.L1].P0];
                TIN_Point P1 = points[lines[A.L1].P1];
                if (TOOL.Existence(P0.Z, P1.Z, h))
                {
                    m = 1;
                    return A.L1;
                }
            }
            if (edgeline != A.L2)
            {
                TIN_Point P0 = points[lines[A.L2].P0];
                TIN_Point P1 = points[lines[A.L2].P1];
                if (TOOL.Existence(P0.Z, P1.Z, h))
                {
                    m = 2;
                    return A.L2;
                }
            }
            return -1;
        }

        //生成等高线
        public void create_contour(double distance)
        {
            if (lines.Count < 3 || points.Count < 3)
            {
                return;
            }
            //获取最小高程和最大高程，代表起始高程与终止高程
            double z_min = points[0].Z, z_max = points[0].Z;
            for (int i = 0; i < points.Count; i++)
            {
                double h = points[i].Z;
                z_min = Math.Min(z_min, h);
                z_max = Math.Max(z_max, h);
            }
            //对所有点的高程进行预处理
            for (double h = Math.Floor(z_min); h < z_max; h += distance)
            {
                for (int i = 0; i < points.Count; i++)
                {
                    if (points[i].Z == h)
                    {
                        points[i].Z += 0.0001;
                    }
                }
            }
            //可用三角形集
            List<int> usabletris = new List<int>();
            #region 处理所有等值线
            //循环所有等值线
            for (double h = Math.Floor(z_min); h < z_max; h += distance)
            {
                usabletris.Clear();
                //生成可用三角形集
                for (int i = 0; i < tris.Count; i++)
                {
                    usabletris.Add(i);
                }
                #region 处理所有边界三角形,生成开曲线
                //处理由边界三角形生成等高线
                for (int i = 0; i < usabletris.Count; i++)
                {
                    if (tris[usabletris[i]].Edge == true)
                    {
                        //获取边界边的编号
                        int edgeline = TOOL.EdgeLine(tris[usabletris[i]]);
                        //判断边界边是否存在等值点
                        TIN_Point P0 = points[lines[edgeline].P0];
                        TIN_Point P1 = points[lines[edgeline].P1];
                        //判断TIN线段上是否存在高度为h的点
                        if (TOOL.Existence(P0.Z, P1.Z, h))
                        {
                            //处理边界三角形，作为等高线的起始
                            Contour contour = new Contour();
                            contour.h = h;
                            //得到边界边上的等值点，并存入该条等值线中
                            TIN_Point temp = TOOL.CalIsoVal(P0, P1, h);
                            contour.points.Add(temp);       
                     
                            int location = -1;//存储另一条等值点所在边位于L0\L1\L2的位置，有0、1、2三种可能，方便找邻接三角形
                            int index = lineindex(edgeline, tris[usabletris[i]],h, ref location);//index为另一个等值点所在边
                            int triindex = tris[usabletris[i]].Adjoin[location];//记录另一条等值点所在边的邻接三角形编号
                            P0 = points[lines[index].P0];
                            P1 = points[lines[index].P1];
                            temp = new TIN_Point();
                            temp = TOOL.CalIsoVal(P0, P1, h);
                            contour.points.Add(temp);
                            usabletris.Remove(usabletris[i]);
                            //追踪等值线
                            while (true)
                            {
                                //如果某条边的相邻三角形的下标为-1，退出循环
                                if (triindex == -1)
                                {
                                    contours.Add(contour);
                                    break;
                                } 
                                int t = triindex;
                                index = lineindex(index, tris[t], h, ref location);
                                triindex = tris[t].Adjoin[location];//记录另一个等值点所在边的邻接三角形编号
                                //计算另一个等值点的坐标，并存入该条等值线中
                                P0 = points[lines[index].P0];
                                P1 = points[lines[index].P1];
                                temp = new TIN_Point();
                                temp = TOOL.CalIsoVal(P0, P1, h);
                                contour.points.Add(temp);
                                //将此三角形从可用三角形数据集中删除
                                usabletris.Remove(t);
                            }
                        }
                        else
                        {
                            usabletris.Remove(usabletris[i]);
                        }
                        i = 0;
                    }
                }
                #endregion

                #region 处理所有闭合等值线
                for (int i = 0; i < usabletris.Count; i++)
                {
                    //获取该三角形的三条边下标号
                    int L0 = tris[usabletris[i]].L0;
                    int L1 = tris[usabletris[i]].L1;
                    int L2 = tris[usabletris[i]].L2;
                    //判断该三角形是否存在等值点，若不存在，将该三角形从可用三角形集中删除
                    if (TOOL.Existence(points[lines[L0].P0].Z, points[lines[L0].P1].Z, h))
                    {
                        //记录该三角形编号
                        int startindex = usabletris[i];
                        //处理该三角形，作为等高线的起始
                        Contour contour = new Contour();
                        contour.h = h;
                        //得到L0边上的等值点，并存入该条等值线中
                        TIN_Point temp = TOOL.CalIsoVal(points[lines[L0].P0], points[lines[L0].P1], h);
                        contour.points.Add(temp);

                        int location = -1;//存储另一个等值点所在边位于L0\L1\L2的位置，有0、1、2三种可能，方便找邻接三角形
                        int index = lineindex(L0, tris[startindex], h, ref location);//index为另一个等值点所在边
                        int triindex = tris[usabletris[i]].Adjoin[location];//记录另一个等值点所在边的邻接三角形编号

                        TIN_Point P0 = points[lines[index].P0];
                        TIN_Point P1 = points[lines[index].P1];
                        temp = new TIN_Point();//temp存储临时等值点
                        temp = TOOL.CalIsoVal(P0, P1, h);
                        contour.points.Add(temp);
                        usabletris.Remove(usabletris[i]);
                        //追踪等值线
                        while (true)
                        {
                            //如果某条边的相邻三角形的下标等于起始三角形的下标，退出循环
                            if (triindex == startindex)
                            {
                                contours.Add(contour);
                                break;
                            }
                            int t = triindex;
                            int tt = index;
                            index = lineindex(tt, tris[t], h, ref location);
                            //index = lineindex(index, tris[t], h, ref location);
                            triindex = tris[t].Adjoin[location];//记录另一个等值点所在边的邻接三角形编号
                            //计算另一个等值点的坐标，并存入该条等值线中
                            P0 = points[lines[index].P0];
                            P1 = points[lines[index].P1];
                            temp = new TIN_Point();
                            temp = TOOL.CalIsoVal(P0, P1, h);
                            contour.points.Add(temp);
                            //将此三角形从可用三角形数据集中删除
                            usabletris.Remove(t);
                        }
                        i = 0;
                    }
                    else if (TOOL.Existence(points[lines[L1].P0].Z, points[lines[L1].P1].Z, h))
                    {
                        int startindex = usabletris[i];
                        //处理该三角形，作为等高线的起始
                        Contour contour = new Contour();
                        contour.h = h;
                        //得到L0边上的等值点，并存入该条等值线中
                        TIN_Point temp = TOOL.CalIsoVal(points[lines[L1].P0], points[lines[L1].P1], h);
                        contour.points.Add(temp);

                        int location = -1;//存储另一个等值点所在边位于L0\L1\L2的位置，有0、1、2三种可能，方便找邻接三角形
                        int index = lineindex(L1, tris[startindex], h, ref location);//index为另一个等值点所在边
                        int triindex = tris[usabletris[i]].Adjoin[location];//记录另一个等值点所在边的邻接三角形编号

                        TIN_Point P0 = points[lines[index].P0];
                        TIN_Point P1 = points[lines[index].P1];
                        temp = new TIN_Point();//temp存储临时等值点
                        temp = TOOL.CalIsoVal(P0, P1, h);
                        contour.points.Add(temp);
                        usabletris.Remove(usabletris[i]);
                        //追踪等值线
                        while (true)
                        {
                            //如果某条边的相邻三角形的下标等于起始三角形的下标，退出循环
                            if (triindex == startindex)
                            {
                                contours.Add(contour);
                                break;
                            }
                            int t = triindex;
                            int tt = index;
                            index = lineindex(tt, tris[t], h, ref location);
                            //index = lineindex(index, tris[t], h, ref location);
                            triindex = tris[t].Adjoin[location];//记录另一个等值点所在边的邻接三角形编号
                            //计算另一个等值点的坐标，并存入该条等值线中
                            P0 = points[lines[index].P0];
                            P1 = points[lines[index].P1];
                            temp = new TIN_Point();
                            temp = TOOL.CalIsoVal(P0, P1, h);
                            contour.points.Add(temp);
                            //将此三角形从可用三角形数据集中删除
                            usabletris.Remove(t);
                        }
                        i = 0;
                    }
                    else
                    {
                        usabletris.Remove(usabletris[i]);
                        i = 0;
                    }
                }
                #endregion
            }
            #endregion
        }

        //计算土方量
        public double[] create_earthwork(double height)
        {
            double[] value = {0,0};//value[0]代表挖方量，value[1]代表填方量
            for (int i = 0; i < tris.Count; i++)
            {
                /*获取第i个三角形，不是L0边上的端点*/
                int test_point = -1;
                if (lines[tris[i].L0].P0 == lines[tris[i].L1].P0 ||
                    lines[tris[i].L0].P1 == lines[tris[i].L1].P0)
                {
                    test_point = lines[tris[i].L1].P1;
                }
                if (lines[tris[i].L0].P0 == lines[tris[i].L1].P1 ||
                    lines[tris[i].L0].P1 == lines[tris[i].L1].P1)
                {
                    test_point = lines[tris[i].L1].P0;
                }
                //获取三角形的三个顶点
                TIN_Point A = points[lines[tris[i].L0].P0];
                TIN_Point B = points[lines[tris[i].L0].P1];
                TIN_Point C = points[test_point];
                //获取三角形与平整面的关系
                int sign = TOOL.Relation(A.Z, B.Z, C.Z, height);
                if (sign == 1)
                {
                    //平整面位于三角形下方
                    value[0] += TOOL.Area(A, B, C) * ((A.Z + B.Z + C.Z) / 3 - height);
                }
                else if (sign == -1)
                {
                    //平整面位于三角形上方
                    value[1] += TOOL.Area(A, B, C) * ((A.Z + B.Z + C.Z) / 3 - height);
                }
                else
                {
                    //平整面位于三角形中间            
                    if (TOOL.Existence(A.Z, B.Z, height))
                    {
                        //线性插值得到平整面与三角形的两个交点,并判断三角形没有等值点的那条边与设计高关系
                        TIN_Point E = TOOL.CalIsoVal(A, B, height);
                        TIN_Point F = new TIN_Point();
                        if (TOOL.Existence(A.Z, C.Z, height))
                        {
                            F = TOOL.CalIsoVal(A, C, height);
                            sign = TOOL.Relation2(B.Z, C.Z, height);
                            double t1 = TOOL.Area(A, E, F) * (A.Z - height) / 3;//该三角形中的被截取的三棱锥部分
                            double t2 = TOOL.Area(A, B, C) * (A.Z - (B.Z + C.Z) / 2) / 3 - t1;//该三角形中的截取后剩下的的棱台部分
                            if (sign == 1)
                            {
                                //说明t1为填方量，t2为挖方量
                                value[0] += (-t2);
                                value[1] += t1;
                            }
                            else
                            {
                                //说明t1为挖方量，t2为填方量
                                value[0] += t1;
                                value[1] += (-t2);
                            }
                        }
                        else
                        {
                            F = TOOL.CalIsoVal(B, C, height);
                            sign = TOOL.Relation2(A.Z, C.Z, height);
                            double t1 = TOOL.Area(B, E, F) * (B.Z - height) / 3;//该三角形中的被截取的三棱锥部分
                            double t2 = TOOL.Area(A, B, C) * (B.Z - (A.Z + C.Z) / 2) / 3 - t1;//该三角形中的截取后剩下的的棱台部分
                            if (sign == 1)
                            {
                                //说明t1为填方量，t2为挖方量
                                value[0] += (-t2);
                                value[1] += t1;
                            }
                            else
                            {
                                //说明t1为挖方量，t2为填方量
                                value[0] += t1;
                                value[1] += (-t2);
                            }
                        }

                    }
                    else if (TOOL.Existence(A.Z, C.Z, height))
                    {
                        //线性插值得到平整面与三角形的两个交点,并判断三角形没有等值点的那条边与设计高关系
                        TIN_Point E = TOOL.CalIsoVal(A, C, height);
                        TIN_Point F = TOOL.CalIsoVal(B, C, height);
                        double t1 = TOOL.Area(C, E, F) * (C.Z - height) / 3;//该三角形中的被截取的三棱锥部分
                        double t2 = TOOL.Area(A, B, C) * (C.Z - (A.Z + B.Z) / 2) / 3 - t1;//该三角形中的截取后剩下的的棱台部分
                        sign = TOOL.Relation2(A.Z, B.Z, height);
                        if (sign == 1)
                        {
                            //说明t1为填方量，t2为挖方量
                            value[0] += (-t2);
                            value[1] += t1;
                        }
                        else
                        {
                            //说明t1为挖方量，t2为填方量
                            value[0] += t1;
                            value[1] += (-t2);
                        }
                    }
                }
            }
            return value;
        }

        //绘制高程点
        public Bitmap Draw_Points(int width, int height,Color color)
        {
            if (points.Count <= 3)
            {
                return new Bitmap(width, height);
            }
            Bitmap image = new Bitmap(width, height);
            Graphics g = Graphics.FromImage(image);
            float min_x = (float)points[0].X, min_y = (float)points[0].Y;
            float max_x = (float)points[0].X, max_y = (float)points[0].Y;
            for (int i = 0; i < points.Count; i++)
            {
                min_x = Math.Min(min_x, (float)points[i].X);
                max_x = Math.Max(max_x, (float)points[i].X);
                min_y = Math.Min(min_y, (float)points[i].Y);
                max_y = Math.Max(max_y, (float)points[i].Y);
            }
            g.TranslateTransform(0, height);
            g.ScaleTransform(1, -1);
            float scale = Math.Max((max_x - min_x) / width, (max_y - min_y) / height);
            g.TranslateTransform(-min_x / scale, -min_y / scale);
            Brush brush = new SolidBrush(color);
            for (int i = 0; i < points.Count; i++)
            {
                float x = (float)points[i].X;
                float y = (float)points[i].Y;
                g.FillEllipse(brush, x / scale - 2, y / scale - 2, 4, 4);
            }
            return image;
        }

        //绘制三角网
        public Bitmap Draw_Lines(int width, int height,Color color)
        {
            if (tris.Count == 0)
            {
                return new Bitmap(width, height);
            }
            Bitmap image = new Bitmap(width, height);
            Graphics g = Graphics.FromImage(image);
            float min_x = (float)points[0].X, min_y = (float)points[0].Y;
            float max_x = (float)points[0].X, max_y = (float)points[0].Y;
            for (int i = 0; i < points.Count; i++)
            {
                min_x = Math.Min(min_x, (float)points[i].X);
                max_x = Math.Max(max_x, (float)points[i].X);
                min_y = Math.Min(min_y, (float)points[i].Y);
                max_y = Math.Max(max_y, (float)points[i].Y);
            }
            g.TranslateTransform(0, height);
            g.ScaleTransform(1, -1);
            float scale = Math.Max((max_x - min_x) / width, (max_y - min_y) / height);
            g.TranslateTransform(-min_x / scale, -min_y / scale);
            Pen pen = new Pen(color, 1);
            for (int i = 0; i < lines.Count; i++)
            {
                PointF p1 = new PointF((float)points[lines[i].P0].X / scale, (float)points[lines[i].P0].Y / scale);
                PointF p2 = new PointF((float)points[lines[i].P1].X / scale, (float)points[lines[i].P1].Y / scale);
                g.DrawLine(pen, p1, p2);
            }
            return image;
        }

        //绘制等高线
        public Bitmap Draw_Contour(int width, int height, Color color)
        {
            if (contours.Count == 0)
            {
                return new Bitmap(width, height);
            }
            Bitmap image = new Bitmap(width, height);
            Graphics g = Graphics.FromImage(image);
            float min_x = (float)points[0].X, min_y = (float)points[0].Y;
            float max_x = (float)points[0].X, max_y = (float)points[0].Y;
            for (int i = 0; i < points.Count; i++)
            {
                min_x = Math.Min(min_x, (float)points[i].X);
                max_x = Math.Max(max_x, (float)points[i].X);
                min_y = Math.Min(min_y, (float)points[i].Y);
                max_y = Math.Max(max_y, (float)points[i].Y);
            }
            g.TranslateTransform(0, height);
            g.ScaleTransform(1, -1);
            float scale = Math.Max((max_x - min_x) / width, (max_y - min_y) / height);
            g.TranslateTransform(-min_x / scale, -min_y / scale);
            Pen pen = new Pen(color, 1);
            PointF[] temp;
            for (int i = 0; i < contours.Count; i++)
            {
                temp = new PointF[contours[i].points.Count];
                for(int j=0;j<contours[i].points.Count;j++)
                {
                    PointF p = new PointF((float)contours[i].points[j].X / scale, (float)contours[i].points[j].Y / scale);
                    temp[j] = p;
                }
                g.DrawLines(pen, temp);
            }
            return image;
        }

        //清空对象数据
        public void Clear()
        {
            points.Clear();
            lines.Clear();
            tris.Clear();
            contours.Clear();
        }

        //清空等高线数据
        public void Clear_Contour()
        {
            contours.Clear();
        }

        //导出TIN为DXF
        public void Export_TIN_DXF(string path)
        {
            StreamWriter sw = new StreamWriter(path);
            sw.WriteLine("0");  // 0在lisp组码中表示开始组码
            sw.WriteLine("SECTION");  // 与后面的ENDSEC对应，类似于{}
            sw.WriteLine("2");
            sw.WriteLine("HEADER");
            sw.WriteLine("0");
            sw.WriteLine("ENDSEC");
            for (int i = 0; i < lines.Count; i++)
            {
                sw.WriteLine("0");
                sw.WriteLine("SECTION");
                sw.WriteLine("2");
                sw.WriteLine("ENTITIES");
                sw.WriteLine("0");
                sw.WriteLine("LINE");
                sw.WriteLine("100");
                sw.WriteLine("AcDbEntity");
                sw.WriteLine("8");
                sw.WriteLine("0");
                sw.WriteLine("100");
                sw.WriteLine("AcDbLine");

                sw.WriteLine("10");  // 起点X坐标
                sw.WriteLine(points[lines[i].P0].X);
                sw.WriteLine("20");  // 起点Y坐标
                sw.WriteLine(points[lines[i].P0].Y);
                sw.WriteLine("30");  // 起点Z坐标
                sw.WriteLine(points[lines[i].P0].Z);
                sw.WriteLine("11");  // 终点X坐标
                sw.WriteLine(points[lines[i].P1].X);
                sw.WriteLine("21");  // 终点Y坐标
                sw.WriteLine(points[lines[i].P1].Y);
                sw.WriteLine("31");  // 终点Z坐标
                sw.WriteLine(points[lines[i].P1].Z);
                sw.WriteLine("0");
                sw.WriteLine("ENDSEC");  // 对应SECTION
            }

            sw.WriteLine("0");
            sw.WriteLine("EOF");  // 文件结束
            sw.Close();  // 关闭文件流
        }

        //导出Point为DXF
        public void Export_Point_DXF(string path)
        {
            StreamWriter sw = new StreamWriter(path);
            sw.WriteLine("0");  // 0在lisp组码中表示开始组码
            sw.WriteLine("SECTION");  // 与后面的ENDSEC对应，类似于{}
            sw.WriteLine("2");
            sw.WriteLine("HEADER");
            sw.WriteLine("0");
            sw.WriteLine("ENDSEC");
            for (int i = 0; i < points.Count; i++)
            {
                sw.WriteLine("0");
                sw.WriteLine("SECTION");
                sw.WriteLine("2");
                sw.WriteLine("ENTITIES");
                sw.WriteLine("0");
                sw.WriteLine("POINT");
                sw.WriteLine("100");
                sw.WriteLine("AcDbEntity");
                sw.WriteLine("8");
                sw.WriteLine("0");
                sw.WriteLine("100");
                sw.WriteLine("AcDbLine");

                sw.WriteLine("10");  // X坐标
                sw.WriteLine(points[i].X);
                sw.WriteLine("20");  // Y坐标
                sw.WriteLine(points[i].Y);
                sw.WriteLine("30");  // Z坐标
                sw.WriteLine(points[i].Z);
                sw.WriteLine("0");
                sw.WriteLine("ENDSEC");  // 对应SECTION
            }

            sw.WriteLine("0");
            sw.WriteLine("EOF");  // 文件结束
            sw.Close();  // 关闭文件流
        }

        //导出等高线为DXF
        public void Export_Contour_DXF(string path)
        {
            StreamWriter sw = new StreamWriter(path);
            sw.WriteLine("0");  // 0在lisp组码中表示开始组码
            sw.WriteLine("SECTION");  // 与后面的ENDSEC对应，类似于{}
            sw.WriteLine("2");
            sw.WriteLine("HEADER");
            sw.WriteLine("0");
            sw.WriteLine("ENDSEC");
            for (int i = 0; i < contours.Count; i++)
            {
                for (int j = 0; j < contours[i].points.Count-1; j++)
                {
                    sw.WriteLine("0");
                    sw.WriteLine("SECTION");
                    sw.WriteLine("2");
                    sw.WriteLine("ENTITIES");
                    sw.WriteLine("0");
                    sw.WriteLine("LINE");
                    sw.WriteLine("100");
                    sw.WriteLine("AcDbEntity");
                    sw.WriteLine("8");
                    sw.WriteLine("0");
                    sw.WriteLine("100");
                    sw.WriteLine("AcDbLine");
                    sw.WriteLine("10");  // 起点X坐标
                    sw.WriteLine(contours[i].points[j].X);
                    sw.WriteLine("20");  // 起点Y坐标
                    sw.WriteLine(contours[i].points[j].Y);
                    sw.WriteLine("30");  // 起点Z坐标
                    sw.WriteLine(contours[i].points[j].Z);
                    sw.WriteLine("11");  // 终点X坐标
                    sw.WriteLine(contours[i].points[j + 1].X);
                    sw.WriteLine("21");  // 终点Y坐标
                    sw.WriteLine(contours[i].points[j + 1].Y);
                    sw.WriteLine("31");  // 终点Z坐标
                    sw.WriteLine(contours[i].points[j + 1].Z);
                    sw.WriteLine("0");
                    sw.WriteLine("ENDSEC");  // 对应SECTION
                }
            }
            sw.WriteLine("0");
            sw.WriteLine("EOF");  // 文件结束
            sw.Close();  // 关闭文件流
        }

        //导出拓扑数据
        public void Export_Topology(string path)
        {
            //输出点
            StreamWriter out_point = new StreamWriter(path+"\\Point.csv",false,Encoding.Default);
            out_point.WriteLine("序号,X坐标,Y坐标,Z坐标");
            for (int i = 0; i < points.Count; i++)
            {
                string data = i.ToString() +","+ points[i].X.ToString("0.000") + ',' ;
                data += points[i].Y.ToString("0.000") + ',' + points[i].Z.ToString("0.000");
                out_point.WriteLine(data);
            }
            out_point.Close();
            //输出线段
            StreamWriter out_line = new StreamWriter(path + "\\Line.csv", false, Encoding.Default);
            out_line.WriteLine("序号,起始点标号,终止点标号");
            for (int i = 0; i < lines.Count; i++)
            {
                string data = i.ToString() + "," + lines[i].P0.ToString() + "," + lines[i].P1.ToString();
                out_line.WriteLine(data);
            }
            out_line.Close();
            //输出三角形
            StreamWriter out_tri = new StreamWriter(path + "\\Tri.csv", false, Encoding.Default);
            out_tri.WriteLine("序号,线段1标号,线段2标号,线段3标号,线段1邻接三角形,线段2邻接三角形,线段3邻接三角形");
            for (int i = 0; i < tris.Count; i++)
            {
                string data = i.ToString() + "," + tris[i].L0.ToString() + "," + tris[i].L1.ToString();
                data += "," + tris[i].L2.ToString() + "," + tris[i].Adjoin[0].ToString() + ",";
                data += tris[i].Adjoin[1].ToString() + "," + tris[i].Adjoin[2].ToString();
                out_tri.WriteLine(data);
            }
            out_tri.Close();
        }

        //获取Point点的个数
        public int Count_Point()
        {
            return points.Count;
        }

        //获取Tri三角形的个数
        public int Count_Tri()
        {
            return tris.Count;
        }

        //获取等高线的条数
        public int Count_Contour()
        {
            return contours.Count;
        }

        //返回离散点数据
        public TIN_Point[] Get_Point()
        {
            TIN_Point[] temp = points.ToArray();
            return temp;
        }

        //返回基线数据
        public TIN_Line[] Get_Line()
        {
            TIN_Line[] temp = lines.ToArray();
            return temp;
        }

        //返回三角形数据
        public TIN_Tri[] Get_Tri()
        {
            TIN_Tri[] temp = tris.ToArray();
            return temp;
        }
    }
}
